#include "vtkMeshImport.hpp"
#include "vtkFunctions.hpp"
#include "vtkFunctions_import.hpp"
#include "vtkConstants.hpp"
#include "meshFunctions.hpp"
#include "meshImportExport.hpp"
#include "meshTypes.hpp"
#include "streamAccess.hpp"
#include "exceptions.hpp"
#include "signalCodes.hpp"

namespace hms::vtk {

namespace xml = tinyxml2;


std::pair<Index, Index> parseExtent(
	const xml::XMLElement* el,
	const char* type
){
	Index nCellsX {0}, nCellsY {0};
	if ( const char* value = el->Attribute( "WholeExtent" ) ){
		std::stringstream converter;
		converter << value;
		Index ignore;
		converter >> ignore >> nCellsX >> ignore >> nCellsY;
		if ( nCellsX == 0 || nCellsY == 0 )
			throw hms::ParseError(
				"(Mesh) Number of cells in x- and y-direction "
				"must not be zero!"
			);
	} else
		throw hms::ParseError(
			"(Mesh) Could not find attribute \"WholeExtent\" in node \""s
			+ type + "\"!"
		);
	return {nCellsX, nCellsY};
}


Matrix2Xs textToVertices( const xml::XMLElement* dataArray, Index cols ){
	/* we're not using read(stringstream, Eigen) on the whole vertex array here,
	 * because VTK always uses 3D coords for vertices, but we only need 2D,
	 * so we're always ignoring every third value.
	 * Alternatively, we could read all values into a Matrix3Xs and then copy
	 * the top rows into a Matrix2Xs. */
	Matrix2Xs vertices (2, cols );
	if (dataArray){
		std::stringstream converter;
		converter << dataArray->GetText();
		scalar ignore;
		for ( Matrix2Xs::ColXpr col : vertices.colwise() ){
			read(converter, col);
			converter >> ignore; /* z coordinate */
		}
	} else
		throw hms::ParseError( "(Mesh) Could not parse vertices!" );
	return vertices;
}

MeshVariant extensionToMeshType( const fs::path& meshPath ){

	MeshVariant meshVar;
	if ( meshPath.extension() == uniExt )
		meshVar.emplace<UniMesh>();
	else if ( meshPath.extension() == rectExt )
		meshVar.emplace<RectMesh>();
	else if ( meshPath.extension() == structExt )
		meshVar.emplace<StructMesh>();
	else if ( meshPath.extension() == unstrExt )
		meshVar.emplace<Mesh>();
	
	return meshVar;
}

void toMesh(
	const xml::XMLDocument* domainDoc,
	const xml::XMLDocument* boundaryDoc,
	MeshVariant& meshVar
){
	toMesh( domainDoc, meshVar);

	readBoundaryPatches( boundaryDoc, meshVar );
}



void toMesh(
	const xml::XMLDocument* doc,
	MeshVariant& meshVar
){
	const xml::XMLElement* vtk { getVtkNode(doc) };
	
	/* call the appropriate function depending on the type held by meshVar */
	meshVar = std::visit(
		[&]( auto&& mesh ) -> MeshVariant {
			using MeshType = remove_qualifiers<decltype(mesh)>;
			const xml::XMLElement* meshNode { getMeshNode<MeshType>(vtk) };
			return toMesh<MeshType>( meshNode );
		},
		meshVar
	);
}


template<>
Mesh  toMesh<Mesh>( const xml::XMLElement* data ){

	/* cells in VTK are represented by their vertices,
	 * so the mesh is reconstructed from them.
	 * First, the vertex coordinates are read. */
	Index nDomainCells {0};
	Matrix2Xs vertices;
	const xml::XMLElement* piece { data->FirstChildElement("Piece") };
	if (piece){
		piece->QueryInt64Attribute(
			( "NumberOf"s + ( usePolyData ? "Polys"s : "Cells"s ) ).c_str(),
			&nDomainCells
		);
		if ( nDomainCells == 0 )
			throw hms::ParseError(
				"(Mesh) NumberOf" + ( usePolyData ? "Polys"s : "Cells"s ) +
				" must not be zero!"
			);
		Index nVertices {0};
		piece->QueryInt64Attribute("NumberOfPoints", &nVertices);
		if (nVertices == 0)
			throw hms::ParseError( "(Mesh) NumberOfPoints must not be zero!" );
		
		const xml::XMLElement* points { piece->FirstChildElement("Points") };
		if (!points)
			throw hms::ParseError(
				"(Mesh) Could not find node \"Points\" in \"Piece\"!"
			);
		if ( const xml::XMLElement* data = points->FirstChildElement("DataArray") )
			vertices = { textToVertices( data, nVertices) };
		else
			throw hms::ParseError(
				"(Mesh) Could not find node \"DataArray\" in \"Points\"!"
			);
	} else 
		throw hms::ParseError(
			"(Mesh) Could not find node \"Piece\" in \""s + unstrType + "\"!"
		);
	
	/* cells are represented as a loop of edges between vertices,
	 * "connectivity" holds the vertex indices which, when drawing lines
	 * between them, form the cells.
	 * "offsets" tell how many vertices form a loop.
	 * For example: Connectivity = 0 1 2 0 2 3, Offsets = 3 6
	 * That means, the first three (=3-0) vertices (0, 1, 2) form a loop (=cell)
	 * as well as the following three (=6-3) */
	const xml::XMLElement* cells { piece->FirstChildElement(
		( usePolyData ? "Polys"s : "Cells"s ).c_str()
	) };

	/* we know the number of offsets, because it's one for each cell.
	 * But we don't know the number of connectivity entries, because there is
	 * more than one way to connect a specific number of vertices
	 * to a specific number of cells */
	ArrayXi connectivity, offsets (nDomainCells);

	if (cells){
		const xml::XMLElement
			* dataBeg { cells->FirstChildElement("DataArray") },
			* data;
		data = dataBeg;
		while (data){
			if ( data->Attribute("Name", "offsets") ){
				textToMatrix( data, offsets);
				break;
			}
			data = data->NextSiblingElement("DataArray");
		}
		/* the last entry in "offsets" equals the number of connectivity entries */
		connectivity.resize( offsets.tail(1)[0] );
		data = dataBeg;
		while (data){
			if ( data->Attribute("Name", "connectivity") ){
				textToMatrix( data, connectivity);
				break;
			}
			data = data->NextSiblingElement("DataArray");
		}

	} else
		throw hms::ParseError(
			"(Mesh) Could not find node \""s +
			( usePolyData ? "Polys"s : "Cells"s ) +
			"\" in \"Piece\"!"
		);

	if ( connectivity.size() == 0 || offsets.size() == 0 )
		throw hms::ParseError(
			"(Mesh) connectivity or offsets not found or empty!"
		);

	/* cells must have either 3 or 4 edges and vertices */
	Array4Xi domainCellVertices (4, nDomainCells );
	ArrayX<char> nDomainCellEdges ( nDomainCells );

	for (int i=0; i<nDomainCells; ++i){
		Index prevOffset { (i==0 ? 0 : offsets[i-1]) };
		Index nEdges = offsets[i] - prevOffset;
		assert( nEdges <= 4 );
		domainCellVertices.col(i).head(nEdges) = connectivity.segment(
			prevOffset, nEdges
		);
		if ( nEdges < 4 )
			domainCellVertices.col(i).tail(4-nEdges) = signalCode::uninitialised;
		nDomainCellEdges[i] = static_cast<char>(nEdges);
	}

	return {
		std::move(vertices), std::move(nDomainCellEdges), domainCellVertices
	};
}



template<>
UniMesh  toMesh<UniMesh>( const xml::XMLElement* data ){

	auto [nCellsX, nCellsY] = parseExtent(data, uniType);

	Vector2s origin;
	origin.array() = std::numeric_limits<scalar>::min();
	if ( const char* value = data->Attribute( "Origin" ) ){
		std::stringstream converter;
		converter << value;
		read(converter, origin);
		if ( (origin.array() == std::numeric_limits<scalar>::min()).any() )
			throw hms::ParseError( "(UniMesh) Could not read origin!");
	} else
		throw hms::ParseError(
			"(UniMesh) Could not find attribute \"Origin\" in node \""s
			+ uniType + "\"!"
		);
	
	scalar dx {std::numeric_limits<scalar>::min()};
	scalar dy {std::numeric_limits<scalar>::min()};
	if ( const char* value = data->Attribute( "Spacing" ) ){
		std::stringstream converter;
		converter << value;
		converter >> dx >> dy;
		if (
			dx == std::numeric_limits<scalar>::min() ||
			dy == std::numeric_limits<scalar>::min()
		) throw hms::IOStreamError( "(UniMesh) Could not read spacing!");
	} else
		throw hms::ParseError(
			"(UniMesh) Could not find attribute \"Spacing\" in node \""s
			+ uniType + "\"!"
		);
	
	return { nCellsX, nCellsY, dx, dy, origin };
}



template<>
RectMesh  toMesh<RectMesh>( const xml::XMLElement* data ){

	auto [nCellsX, nCellsY] = parseExtent(data, rectType);

	VectorXs xCoords (nCellsX+1);
	VectorXs yCoords (nCellsY+1);
	if (
		const xml::XMLElement* coords =
		data->FirstChildElement("Piece")->FirstChildElement("Coordinates")
	){
		textToMatrix( coords->FirstChildElement("DataArray"), xCoords);
		textToMatrix(
			coords->FirstChildElement("DataArray")->
				NextSiblingElement("DataArray"),
			yCoords
		);
	} else
		throw hms::ParseError(
			"(RectMesh) Could not find node \"Coordinates\" in node \""s
			+ rectType + "\"->\"Piece\"!"
		);
	
	return { std::move(xCoords), std::move(yCoords) };
}



template<>
StructMesh  toMesh<StructMesh>( const xml::XMLElement* data ){

	auto [nCellsX, nCellsY] = parseExtent(data, structType);

	const xml::XMLElement* dataArray {
		getChildNodeByName(
			getChildNodeByName(
				getChildNodeByName(data, "Piece"),
				"Points"
			),
			"DataArray"
		) };

	Matrix2Xs vertices { textToVertices(
		dataArray, (nCellsX+1)*(nCellsY+1)
	) };
	
	return { std::move(vertices), nCellsX, nCellsY };
}


void readBoundaryPatches( const xml::XMLNode* parent, MeshVariant& meshVar ){

	const xml::XMLElement* boundary { getBoundaryNode(parent) };

	std::visit( [&]( auto& mesh ){
		using MeshType = remove_qualifiers<decltype(mesh)>;

		Array3Xi edgesAtBoundaryVerts;
		if constexpr ( std::is_same_v<MeshType, Mesh> ){
			edgesAtBoundaryVerts = edgesAtBoundaryVertices(
				mesh.edges().vertices(), mesh.nInnerEdges()
			);
		}
		
		const xml::XMLElement* patch {
			boundary->FirstChildElement(boundaryChild)
		};
		while (patch){
			const char* patchName { patch->Attribute("Name") };
			
			if ( patchName ){
				if ( std::string_view{patchName}.empty() )
					throw hms::ParseError(
						"Boundary patch name must not be empty!"
					);
			} else {
				throw hms::ParseError(
					"Attribute \"name\" in boundary patch node not found!"
				);
			}

			Array3i verts;
			textToMatrix( patch->FirstChildElement("Vertices"), verts );

			addBoundaryPatch( mesh, patchName, verts, edgesAtBoundaryVerts );

			patch = patch->NextSiblingElement(boundaryChild);
		}
	}, meshVar);
}

void readBoundaryPatches( const fs::path& path, MeshVariant& meshVar ){
	std::unique_ptr<xml::XMLDocument> doc { readXML(path) };
	readBoundaryPatches(doc.get(), meshVar);
}

}